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ABSTRACT 


2.  SYSTEM  MODEL  AND  ML  ESTIMATION 


One  main  drawback  of  using  Time  Difference  of  Arrival  (TDOA) 
methods  for  source  localization  and  navigation  is  that  they  require 
centralization  of  multiple  copies  of  a  signal.  This  paper  considers 
blindly  estimating  the  location  of  a  Cyclic  Prefix  (CP)  in  an  Orthog¬ 
onal  Frequency  Division  Multiplexing  (OFDM)  signal,  enabling  dis¬ 
tributed  TDOA  computation  up  to  an  integer  ambiguity.  This  ambi¬ 
guity  can  be  resolved  using  integer  least-squares  methods,  if  enough 
TDOAs  are  available.  The  contributions  of  this  paper  are  derivation 
of  the  Cramer-Rao  Lower  Bound  (CRLB)  on  locating  the  CP,  and 
hence  on  the  underlying  source  localization  problem. 

Index  Terms —  Cramer-Rao  lower  bound,  OFDM,  cyclic  prefix, 
time  delay  estimation,  blind,  TDOA 


1.  INTRODUCTION 

One  main  drawback  of  using  TDOA  methods  for  source  localization 
and  navigation  is  that  they  require  centralization  of  multiple  copies 
of  a  signal,  which  wastes  bandwidth  and  power.  For  OFDM  sources, 
the  amount  of  required  centralization  of  data  can  he  greatly  reduced 
by  comparing  the  temporal  locations  of  the  CPs  rather  than  compar¬ 
ing  the  entire  signals  [1].  This  can  be  done  blindly,  without  knowl¬ 
edge  of  the  data  contained  in  each  CP  [2],  However,  since  CPs  occur 
at  regular  intervals,  this  leaves  an  integer  ambiguity  in  each  TDOA 
(any  integer  times  the  OFDM  block  length).  In  [1],  the  ambigu¬ 
ity  was  resolved  by  transmitting  a  small  amount  of  data  per  block, 
whereas  in  this  paper,  we  assume  that  the  ambiguity  can  be  resolved 
using  integer  least  squares  methods,  as  in  [3], 

Many  papers  discuss  the  CRLB  of  TDOA-based  positioning  [4], 
Existing  work  assumes  cross-correlation  either  between  sensors  [5], 
[6],  [7]  or  with  a  known  training  signal  [8],  In  this  paper,  we  consider 
the  TDOA  being  estimated  blindly  and  in  a  distributed  fashion,  by 
comparing  the  delays  of  the  CPs  (estimated  as  in  [2])  in  pairs  of  re¬ 
ceived  signals.  Thus,  the  main  contribution  of  this  paper  is  derivation 
of  the  CRLB  on  estimating  the  time  delay  of  the  CP  within  a  sym¬ 
bol.  Assuming  that  the  integer  ambiguity  is  appropriately  resolved, 
this  allows  derivation  of  a  CRLB  on  the  TDOAs,  and  ultimately  the 
CRLB  for  the  underlying  source  localization  problem,  for  the  case 
in  which  no  data  is  centralized  aside  from  the  reception  times  of  the 
CPs. 
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The  model  for  the  received  signal  is 

y(t)  =  E  x[k]p(t  —  5  kT)  +  n(t).  (1) 

k 

where  5  is  the  time  delay  to  he  estimated,  x[k\  is  the  transmitted  data 
sequence,  and  p(t)  is  the  pulse  shape  with  derivative  p'(t).  We  will 
use  p0  (f )  to  denote  a  time-normalized  version  of  the  pulse  shape,  i.e. 
p(t)  =  p0(t/T),  and  we  assume  that  p(t)  is  a  raised-cosine  pulse 
with  excess  bandwidth  (3 . 

The  OFDM  transmitter  uses  an  FFT  size  of  N  and  a  CP  length  of 
v  samples.  The  transmitter  sends  one  block  of  M  =  N  +  is  samples 
of  x[k]  every  MT  seconds.  The  receiver  uses  an  oversampling  factor 
of  q ,  and  thus  has  a  sampling  period  of  Ts  =  T/q.  The  receiver 
samples  a  total  of  L  blocks,  or  LqM  samples.  The  set  Tc  denotes 
the  set  of  transmitted  samples  of  x[k]  in  cyclic  prefixes  (the  first  u  of 
each  M  samples),  the  set  Id  denotes  the  unique  N  —  v  samples  in 
the  middle  of  each  M  samples,  and  the  set  Ie  denotes  the  last  v  of 
each  M  samples  (which  get  copied  into  the  CP). 

The  signal  power  and  noise  power  are  and  <r(j.  per  sample, 
respectively,  and  the  Signal  to  Noise  Ratio  (SNR)  is  cr^/cr^.  The 
transmitted  data  is  assumed  to  he  white  (aside  from  the  presence  of 
the  CP),  and  the  noise  is  assumed  to  be  completely  white  over  the 
spectrum  of  interest. 

In  [2],  a  Maximum  Likelihood  (ML)  algorithm  was  derived  to 
jointly  estimate  the  temporal  location  of  the  CP  and  the  Carrier  Fre¬ 
quency  Offset  (CFO),  hut  no  CRLB  was  derived.  In  this  paper,  for 
simplicity,  we  assume  there  is  no  CFO.  The  ML  algorithm  of  [2], 
generalized1  to  allow  for  oversampling  by  q  and  averaging  over  L 
blocks,  is  given  by 
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l{5) 

m 

P 


arg  max  |7(<5)  -  ^  $(5)} ,  (2) 

L  Mql-\-8-\-uq  —  l 

E  E  y(kTs)y*((k  +  Nq)Ts),  (3) 

1=1  k=Mql-\-8 


L  Mql-\-8-\-uq—l 

E  E  {\y(kTs)\2  +  \y((k  +  Nq)Ts)\2),  (4) 

1=1  k=Mql-\-8 


SNR 
SNR + 1 ' 


(5) 


*The  ML  algorithm  of  [2]  assumed  q  =  1,  and  it  is  not  immediately 
obvious  whether  the  generalization  here  is  still  ML  when  q  >  1. 
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including  suggestions  for  reducing  this  burden,  to  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington 
VA  22202-4302.  Respondents  should  be  aware  that  notwithstanding  any  other  provision  of  law,  no  person  shall  be  subject  to  a  penalty  for  failing  to  comply  with  a  collection  of  information  if  it 
does  not  display  a  currently  valid  OMB  control  number. 


1.  REPORT  DATE 

MAY  2010 


2.  REPORT  TYPE 

N/A 


3.  DATES  COVERED 


5a.  CONTRACT  NUMBER 


5b.  GRANT  NUMBER 


5c.  PROGRAM  ELEMENT  NUMBER 


5d.  PROJECT  NUMBER 


5e.  TASK  NUMBER 


5f.  WORK  UNIT  NUMBER 


4.  TITLE  AND  SUBTITLE 

Bounds  On  Distributed  Tdoa-Based  Localization  Of  Of  dm  Sources 


6.  AUTHOR(S) 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES)  8.  PERFORMING  ORGANIZATION 

The  Air  Force  Institute  of  Technology  Dept,  of  Elec.  &  Comp.  Eng.  report  number 

WPAFB,  OH,  45433 

9.  SPONSORING/MONITORING  AGENCY  NAME(S )  AND  ADDRESS(ES )  10.  SPONSOR/MONITOR' S  ACRONYM(S) 

11.  SPONSOR/MONITOR'S  REPORT 
NUMBER(S) 

12.  DISTRIBUTION/AVAILABILITY  STATEMENT 

Approved  for  public  release,  distribution  unlimited 

13.  SUPPLEMENTARY  NOTES 

See  also  ADM002314.  IEEE  International  Conference  on  Acoustics,  Speech  and  Signal  Processing  (34th) 
held  in  Taipei,  Taiwan  on  19-24  April  2009.  U.S.  Government  or  Federal  Purpose  Rights  License,  The 
original  document  contains  color  images. 

14.  ABSTRACT 

One  main  drawback  of  using  Time  Difference  of  Arrival  (TDOA)  methods  for  source  localization  and 
navigation  is  that  they  require  centralization  of  multiple  copies  of  a  signal.  This  paper  considers  blindly 
estimating  the  location  of  a  Cyclic  Prefix  (CP)  in  an  Orthogonal  Frequency  DivisionMultiplexing  (OFDM) 
signal,  enabling  distributed  TDOA  computation  up  to  an  integer  ambiguity.  This  ambiguity  can  be 
resolved  using  integer  least-squares  methods,  if  enough  TDOAs  are  available.  The  contributions  of  this 
paper  are  derivation  of  the  Cramer-Rao  Lower  Bound  (CRLB)  on  locating  the  CP,  and  hence  on  the 
underlying  source  localization  problem. 


15.  SUBJECT  TERMS 


16.  SECURITY  CLASSIFICATION  OF: 


a.  REPORT 

unclassified 


b.  ABSTRACT 

unclassified 


c.  THIS  PAGE 

unclassified 


17.  LIMITATION  OF 

18.  NUMBER 

ABSTRACT 

OF  PAGES 

SAR 

4 

19a.  NAME  OF 
RESPONSIBLE  PERSON 


Standard  Form  298  (Rev.  8-98) 

Prescribed  by  ANSI  Std  Z39-18 


In  the  next  section,  we  will  also  consider  an  approximate  ML  algo¬ 
rithm,  given  by 


Sapx  =  arg  max  {7(6)}  •  (6) 

_M_P<s<M_p 

This  is  motivated  by  the  fact  that  <f>(<5)  is  approximately  constant 
over  <5,  especially  if  the  amount  of  averaging  (L)  is  large.  The  p  <E>(<5) 
term  also  vanishes  at  low  SNR. 

3.  EFFECTS  OF  SAMPLING 


successive  set  of  N  elements  corresponds  to  a  block  of  (unprefixed) 
independent  data  samples  x[k\. 

The  (scalar)  submatrix  A,  evaluated  at  S  =  0  for  simplicity,  is 


\d2C\  1  v  (8E{y(mT)}\ 

1  dS*  j  a2  ^  \  85  ) 


z[m ]  =  x[k]p'((m  —  k)T). 


(11) 

(12) 

(13) 


Intuitively,  the  CRLB  should  depend  on  the  signal  bandwidth  1/T, 
but  not  on  the  actual  sampling  rate  1/TS,  since  the  latter  in  part  de¬ 
fines  the  estimation  method.  However,  the  performance  of  a  given 
estimator  will  depend  on  the  oversampling  factor.  For  example,  the 
ML  algorithm  of  [2]  does  not  interpolate  between  samples,  hence  its 
variance  is  at  least  that  of  a  uniform  random  variable  distributed  over 
[— Ts/2,  Ts/2],  which  is  T2  / 12.  This  “sample  resolution”  bounds 
the  Root  Mean  Squared  Error  (RMSE),  causing  a  departure  from  the 
CRLB  at  high  SNR  values. 

In  order  for  the  CRLB  to  be  analytically  tractable  in  a  short 
paper,  the  derivation  in  the  next  section  assumes  sampling  at  the 
Nyquist  rate,  i.e.  q  =  1,  which  makes  adjacent  samples  uncorre¬ 
lated.  Prior  work  on  time  delay  estimation  [5],  [6]  avoids  this  as¬ 
sumption  by  working  in  the  frequency  domain.  However,  it  does  not 
make  sense  to  discuss  the  power  spectra  of  OFDM  systems,  since  the 
CP  makes  OFDM  non- stationary  (in  fact,  cyclo-stationary  with  pe¬ 
riod  MT).  However,  similar  work  on  image  registration  [9]  makes 
an  assumption  analogous  to  ours,  i.e.  that  the  observed  images  are 
sampled  at  the  Nyquist  rate  or  greater. 

4.  CRLB  DERIVATION 

First,  we  consider  the  general  case,  then  we  derive  the  bound  on  es¬ 
timators  that  only  use  the  autocorrelation  function  of  (3)  to  locate 
the  CP.  The  latter  derivation  is  cleaner  but  only  applies  to  auto¬ 
correlation  based  estimators,  e.g.  (6).  Real  data  is  assumed  for  sim¬ 
plicity,  and  the  SNR  is  assumed  known. 

4.1.  General  CRLB 

In  the  general  approach,  the  unknowns  are  the  time  delay  <5  and  the 
LN  nuisance  parameters  {a:[fc]}.  Let  x  consist  of  L  blocks  of  N 
samples  of  x[k ]  each,  containing  no  repetition,  i.e.  omitting  the  CP 
samples.  The  log-likelihood  of  the  received  vector  y  is 

y  =  ,y(LMT)]T  (7) 

ln(/(y|x,<)))_  =  const.  -  ^  [y(mT)  -  E  {y(mT)}}2  , 

v  71  m 

C 

(8) 

E  {y(mT)}  =  ^  x[k]p(mT  —  5  —  kT).  (9) 

k 

The  Fisher  Information  Matrix  (FIM)  has  a  block  structure, 


Again  for  S  =  0,  B  and  C  are  given  element-wise  by 

^2  r  1  —  8E  {y(mT)}  <9E  {y(mT)} 
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d2C  \  _  J_  V 

cW  8x[ko\  J  a2 


85  8x[ko\ 

8E{y(mTs)} 


-E 


d2C 


{  dx[ko\  dx[mo\ 


8x[k0\ 

8E  {y(mT)}  8E  {y(mT)} 


(14) 

(15) 


dx[k0] 


Due  to  the  repetition  induced  by  the  CP, 
8E{y(mT)} 


8x[k0\ 


f  p((m-  ko)T), 

=  {  p({m  -  ( k0  -  N))T) 
[  +p((m-  k0)T), 


dx\j0\ 


ko  £  Td 


k()  £  IZe 


(16) 


(17) 


For  a  Nyquist  pulse  shape,  (17)  becomes 


8E{y(rnT)} 

8x[ko\ 


5k  (jn,  k0), 


ko  £  Td 


fco  —  N)  +  <5x(m,fco),  ko  €  le 


where  5k(-)  is  the  Kronecker  delta  function.  Thus, 

A= 4i>2h 


B[fc0]  = 


C  = 


—z[ko\,  ko&Td 

—  z[fco  —  N]  —  z[ko\,  ko  e  Te 


In-,  0 

0  2I„ 


I  I L, 


(18) 

(19) 

(20) 
(21) 


where  0  is  the  Kronecker  product.  (Note  that  ko  only  indexes  the 
last  N  of  each  N  +  v  samples  of  x.)  Using  Schur  complements  to 
perform  the  matrix  inversion  in  block  fashion,  the  top  left  element  of 
the  CRLB  is  given  by 


VAR[5]  >  [a-BtC_1b] 


2 

— 


^V[m] 


mexd 


(22) 

-  -1 

l  J2  (ZH  +  z[m  -  N])2 

m£le 

(23) 


J  = 


A 

B 


Bt 
C  ’ 


(10) 


where  A  is  1  x  1,  B  is  LN  x  1,  and  C  is  LN  x  LN.  In  each 
dimension,  the  first  element  of  J  corresponds  to  <5,  and  then  each 


In  (23),  there  are  three  summations  over  the  index  m.  The  first  sum¬ 
mation  includes  all  m,  including  the  CP  set,  Lc;  the  data  set,  Xd\  and 
the  data  in  the  ends  of  blocks,  Te.  The  second  summation  only  in¬ 
cludes  the  data  in  the  middle  of  each  block,  and  the  third  summation 
only  includes  the  data  in  the  ends  of  blocks. 
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Noting  that  “m  £  Ie”  is  equivalent  to  “(m  —  N)  £  lc 


VAR[6 ]  >  2cr„ 


=  2<j2  T2 


J2  {z[m]  -  z[m  -  iV])2 

J2  ~  x^k  ~  N^p'^m  ~ k^] 

nGXe  \  k  ) 


(24) 


(25) 


Observe  that  the  terms  when  k  £  Xe  are  zero.  Moreover,  we  al¬ 
ways  have  m  £  Ie,  and  the  pulse  shape  factor  is  nearly  zero  except 
when  | m  —  k\  is  small.  Thus,  the  terms  that  contribute  most  to  the 
summation  are  the  boundary  terms,  i.e.  when  k  is  just  outside  ofTe. 

Equation  (25)  is  the  CRLB,  and  cannot  be  simplified  further 
without  approximations.  However,  in  order  to  gain  intuition,  con¬ 
sider  the  case  of  large  L,  which  enables  the  approximation 


E(i)  -  <26) 

m£Te  l—l  rn£X*  ^  m^X\  J 


Given  the  triangular  shape  of  r(t),  its  derivative  has  a  magnitude  of 
La2  /T  for  a  range  of  2 u  samples,  and  is  zero  otherwise.  Thus, 


{2v)(Lal/T)2 
Lv(al  +alf 


Inverting  and  simplifying, 


V AR[5Corr\  > 


(l  +  SNR-1)2  T2 
2 L 


The  limiting  cases  are 


VAR[6corr] 


rjp  2 

2 l  SNR2 ' 


high  SNR 
low  SNR 


(32) 


(33) 


(34) 


Interestingly  (and  counter-intuitively),  both  (27)  and  (33)  are  in¬ 
dependent  of  the  fraction  of  each  block  consisting  of  the  CP,  i.e.  the 
fraction  zz/M.  Thus,  a  short  CP  is  as  good  as  a  long  CP  for  purposes 
of  blind  delay  estimation  in  the  absence  of  multipath. 


where  the  set  Ile  is  the  data  in  block  l  of  Ie.  Then,  given  that  the 
data  x[k\  is  uncorrelated  (aside  from  the  repetition  in  the  CP), 


VAR[5]  >  2<t2T2 


l  (2cr^)  - fc))2 

mell  k£Te 


-1 


rp  2 

L  SNR 


rp2 

L  SNR 


OO 

2  ^d(p'o{d))2 

.  d=  1 
v 


(27) 


v 


For  raised  cosine  pulse  shapes  with  excess  bandwidths  /3  of  0, 
0.25,  and  0.5,  the  factor  77  is  0.32,  0.41,  and  0.74,  respectively. 


4.2.  Bound  for  autocorrelation-based  methods 

This  section  derives  the  CRLB  for  estimators  that  only  use  the  auto¬ 
correlation  data  7 (t),  as  in  (6).  Since  7(f)  are  the  observations,  the 
only  unknown  is  the  time  delay  5.  The  expected  value  of  7 (t)  is 

T(f)  =  E  {7(f)}  =  L  a2  A(f  —  S),  (28) 

where  A(t)  is  a  triangular  pulse  spanning  —  v  <t<v  with  a  peak 
height  of  v.  It  can  be  shown  that  the  variance  of  7  (f)  is 


4.3.  CRLB  on  source  localization 


Assume  K  synchronized  TDOA  receivers  work  together  to  locate 
one  OFDM  source.  Given  the  estimates  of  the  CP  locations  and  a 
suitable  method  for  resolving  the  integer  ambiguities  [3],  the  TDOA 
between  each  pair  of  receivers  can  be  obtained  by  subtracting  the 
two  estimated  delays.  This  doubles  the  variance  and  the  CRLB,  and 
induces  some  covariance  between  TDOA  estimates.  Thus,  the  co- 
variance  matrix  of  efficient  TDOA  estimates  and  its  inverse  are 
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Q”1 


77  T2 
L  SNR 


[I  +  111 


( ,T2 1 

-1 

1  t! 

I - 11T 

\l  snr  y 

K  \ 

(35) 

(36) 


where  1  is  a  ( K  —  1)  x  1  vector  of  ones. 

From  [7],  the  CRLB  on  the  position  estimate  can  be  computed 
given  the  covariance  matrix  of  a  TDOA  estimator,  as 


CRLB  (. xtx,ytx )  =  c2 


(gtQ“1g) 


2  rri2 

77  c  T 
L  SNR 


(37) 

(38) 


where  G  is  a  ( K  —  1)  x  2  matrix  that  depends  entirely  on  the 
geometries  of  the  transmitter  and  receivers,  with  the  ith  row  (for 
i  =  2, 3,  ■  ■  ■  ,  K)  given  by 


VAR{1(t)}  =  Lv(al+al)'2 ,  (29) 


and  that  the  auto-covariance  of  7 (ti)  and  7(42)  is  zero  for  all  pairs 
of  non-equal  time  instants. 

Let  7  =  [jy(— M/2),--  -  —  1)].  The  log-likelihood 

function  of  this  vector  given  the  delay  5  is 


ln(/(7|<5))  =  const. 


c 


Et  b(t)  -  mf 

2 Lv  (<ri  +  cr2)2 


The  Fisher  information  is  given  by 


rfa2£)  EtTO] 
l  952  /  Lv  (ct2  +  ■ 


(30) 


(31) 


_  xi  -  xtx  xi  —  xtx  yi  -  ytx  yi  -  ytx 

L*(rowi)  — - 1 - 

77  ri  n  ri 

(39) 

where  (xtx,  ytx )  and  (xt,  yi)  are  the  locations  of  the  source  and  K 
receivers,  and  ri  is  the  distance  from  the  source  to  receiver  i. 

5.  SIMULATIONS 

The  simulation  parameters  in  this  section  are  comparable  to  an  IEEE 
802.11a  system:  the  FFT  size  is  N  =  64,  the  CP  length  is  u  =  16, 
the  block  size  is  M  =  N+  v  =  80,  and  the  Nyquist  sampling  period 
is  T  =  50  ns.  The  number  of  blocks  is  L  =  10  (yielding  a  40  /rs 
observation  time).  The  raised  cosine  pulse  shape  uses  (3  =  0. 
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Fig.  3.  Performance  of  position  estimation,  as  in  Section  4.3,  for 
0  dB  SNR  and  q  =  7.  Each  V  is  a  sensor,  the  *  is  the  source, 
and  each  +  is  a  resolved  source  estimate  (200  trials).  The  outer  and 
inner  ellipses  indicate  the  RMSE  and  the  CRLB,  both  scaled  up  by  2 
so  that  the  former  gives  an  86%  confidence  interval. 


Fig.  1.  Comparison  of  the  CRLBs  of  (27)  and  (33)  to  the  RMSEs  of 
(2)  and  (6).  The  upper  bound  corresponds  to  a  uniform  distribution 
across  a  span  of  MT  seconds,  and  “sample  res”  is  the  floor  discussed 
in  Section  3.  There  is  no  oversampling  ( q  =  1). 


Fig.  2.  Comparison  of  the  CRLB  of  (27)  to  the  RMSE  of  (2)  for 
various  oversampling  factors  q.  The  CRLB  was  derived  assuming 
q  =  1,  but  appears  to  hold  for  all  q. 


Fig.  1  shows  the  square  roots  of  the  CRLBs  from  (27)  and  (33), 
and  the  RMSEs  of  (2)  and  (6),  using  1500  trials.  The  RMSE  floor 
is  as  discussed  in  Section  3.  Fig.  2  compares  the  performance  of  the 
ML  estimator  for  various  values  of  q  to  the  CRLB.  Although  the 
CRLB  was  derived  for  q  =  1,  it  appears  to  hold  for  all  q.  Sampling 
faster  initially  improves  the  performance  of  the  ML  estimator,  but 
the  realm  of  diminishing  returns  is  reached  very  quickly.  Fig.  3  and 
Table  1  compare  the  final  RMSE  and  CRLB  on  position  estimation, 
for  q  =  7.  Ambiguity  resolution  was  performed  by  dividing  the 
space  into  3  km  x  3  km  regions  and  bounding  the  integer  offsets 
within  each  region  based  on  the  sensor  geometry  and  region  size. 
The  best  choice  of  integers  for  each  region  was  determined  using 
box-constrained  integer  least  squares  [3],  The  overall  optimum  inte¬ 
gers  were  determined  by  selecting  each  region’s  choice  in  turn,  esti¬ 
mating  the  position  via  [10],  and  comparing  the  TDOAs  that  would 
result  from  that  estimate  to  the  observations. 


Table  1.  Effects  of  varying  SNR  on  position  estimation  (in  meters). 


SNR: 

-10  dB 

-5  dB 

0  dB 

5  dB 

10  dB 

30  dB 

RMSE 

5790 

11.5 

3.82 

2.26 

1.85 

1.47 

VCRLB 

50.7 

5.07 

1.67 

1.01 

0.81 

0.66 
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